%
% Author: Vasyl Mykhalchuk
%
function [] = scalarFieldFromOBJSequence(averageYES)
format short;

% Modify matlab path
addpath(genpath('../util'));
addpath(genpath('../io'));
addpath(genpath('../strain'));
addpath(genpath('../export'));
addpath(genpath('../obj-tools'));
addpath(genpath('..'));

LOG('Loading dependencies');

%Load OBJ files for reference/target meshes
%{
referencePath = strcat('examples\example3.obj');
targetPath = strcat('examples\example3.obj');
%}
%{
referencePath = strcat('horse\horse-gallop-01.obj');
targetPath = strcat('horse\horse-gallop-16.obj');
%}


referencePath = strcat('..\..\res\head\female_01.obj');
targetPath = strcat('..\..\res\head\female_02.obj');

%{
referencePath = strcat('female4\female4_01.obj');
targetPath = strcat('female4\female4_300.obj');
%}

[vertices, faces] = read_wobj(referencePath);
obj_ref.vertices = vertices;
obj_ref.faces = faces;
[nF,x] = size(faces);
[nV,x] = size(vertices);
obj_ref.nF = nF;
obj_ref.nV = nV;

[vertices, faces] = read_wobj(targetPath);
obj_target.vertices = vertices;
obj_target.faces = faces;

numberOfFaces = obj_ref.nF;
numberOfVertices = obj_ref.nV;

% compute adjacency lists
% List of faces adjacent to each vertex
LOG('Preprocessing');
adjList = createAdjacencyList(obj_ref.faces, numberOfVertices, numberOfFaces);
% List contains list of neighboring vertices for each vertex
adjVertList = createVertexAdjacencyList(obj_ref.faces, numberOfVertices, numberOfFaces);

% time-varying scalar field (vertex based)
scalarFields = zeros(numberOfVertices, 1, 'double');
% time-varying scalar field (triangle based)
scalarFieldsT = zeros(numberOfFaces, 1, 'double');

faces = obj_ref.faces;

% =============== ============== ===========
LOG('Computing faces scalar field');
for iFace = 1:numberOfFaces
    tr = faces(iFace,:);
    v1 = obj_ref.vertices(tr(1), 1:3);
    v2 = obj_ref.vertices(tr(2), 1:3);
    v3 = obj_ref.vertices(tr(3), 1:3);
    
    vv1 = obj_target.vertices(tr(1), 1:3);
    vv2 = obj_target.vertices(tr(2), 1:3);
    vv3 = obj_target.vertices(tr(3), 1:3);
    
    %flatten triangles 
    [v1f, v2f, v3f] = fitTrianglesToOnePlane(v1, v2, v3);
    [vv1f, vv2f, vv3f] = fitTrianglesToOnePlane(vv1, vv2, vv3);
    
    strains = tri2D_strain(v1f, v2f, v3f, vv1f, vv2f, vv3f);
    kTreshold = 0.009;
    if( abs(strains(1)) > kTreshold )
        scalarFieldsT(iFace) = strains(1);
    else 
        scalarFieldsT(iFace)=0;
    end%if
end%for


% compute vertex-based time-varying scalar field
LOG('Computing vertex scalar field');
deformationMax = 0.0;
deformationMin = 0.0;
for iVertex = 1:numberOfVertices
        i = 1;
        scalarValue = 0.0;
        adjTriangle = adjList(iVertex,i);
        while (adjTriangle ~= 0)
            tr = faces(adjTriangle,:);
            v1 = obj_target.vertices(tr(1), 1:3);
            v2 = obj_target.vertices(tr(2), 1:3);
            v3 = obj_target.vertices(tr(3), 1:3);
            trArea = 1; %triangleArea(v1,v2,v3);
            scalarValue = scalarValue + trArea * scalarFieldsT(adjTriangle);
            i = i + 1;
            adjTriangle = adjList(iVertex,i);
        end
        if scalarValue > deformationMax
            deformationMax = scalarValue;
        end%if
        if scalarValue < deformationMin
            deformationMin = scalarValue;
        end%if
        scalarFields(iVertex) = scalarValue;
end%for

LOG('Deformation MAX'); deformationMax
LOG('Deformation MIN'); deformationMin


%LOG('Gathering color info');
%vertexColorList = createColorInfo(linspace(1,numberOfVertices, numberOfVertices), scalarFields);
LOG('Extracting critical points');
cp = extractCriticalPoints(scalarFields, adjVertList);

format long;
scalarFields(cp)


hist(scalarFields(cp));
plot(scalarFields(cp));

for i=1:numberOfVertices
        if deformationMax ~= 0
            if(deformationMax >= abs(deformationMin)) absDeformationMax = deformationMax; else absDeformationMax = abs(deformationMin);end
            %scalarFields(i) = scalarFields(i) / absDeformationMax;
        end
end%for

%LOG('Export scalar field mesh');
%exportPLY( 'target.ply', obj_ref.vertices, vertexColorList, faces);
vertexColorList = createColorInfo(cp, scalarFields);
LOG('Export critical points');
exportPLY( 'target_cp.ply', obj_ref.vertices, vertexColorList, faces);

end%function scalarFieldFromOBJSequence







